{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "author ： duanxxnj@163.com\n",
      "time : 2016-06-06_15-41\n",
      "\n",
      "Lasso 回归应用于稀疏信号\n",
      "\n",
      "\n",
      "测试集上的R2可决系数 : -0.185315\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3Xl8VNX5+PHPmTV7CEnYCYuCrCFAWEQBFQVXUOpOLdSFUtfWn2ttq+3XtrZa963U4r7WXUGLKApuICDIvggIgQCBEELW2c7vjzsz2TdmMneSed6vV153MnfmzjN3Zp45c85zz1Vaa4QQQrQfFrMDEEIIEV6S2IUQop2RxC6EEO2MJHYhhGhnJLELIUQ7I4ldCCHaGUnsQgjRzjQ7sSul5imlDiil1lW7rqNS6hOl1Fb/Mq11whRCCNFcLWmxPwecWeu6O4BPtdb9gE/9/wshhDCRasmRp0qp3sCHWush/v83A6dorfOVUl2Bz7XWJzS1nYyMDN27d+9jClgIIWLVypUrD2qtM5u6nS3Ex+mstc4H8Cf3Ts25U+/evVmxYkWIDy2EELFFKfVTc24XscFTpdRspdQKpdSKgoKCSD2sEELEnFAT+35/Fwz+5YGGbqi1nqu1ztVa52ZmNvlLQgghxDEKNbG/D8z0X54JvBfi9oQQQoSo2X3sSqlXgVOADKVUHnA3cB/whlLqKmAXcFFrBClEe+V2u8nLy6OiosLsUEQUiYuLo0ePHtjt9mO6f7MTu9b6sgZWTTqmRxZCkJeXR3JyMr1790YpZXY4IgporTl06BB5eXn06dPnmLYhR54KYaKKigrS09MlqYsgpRTp6ekh/YqTxC6EySSpi9pCfU9IYheitfl88P1L4HWbHYmIEZLYhWhta16F966Drx8zO5J6JSUlmR1CyG699VYGDx7MrbfeSkFBAWPGjGH48OEsXbqUs88+m6Kiogbv+/TTT/PCCy8c0+Pu3LmTV1555VjDbjWhHnkqhGhKeaGxLDtkbhzt2L/+9S8KCgpwOp289tprDBgwgOeffx6A8ePHN3rfOXPmHPPjBhL75ZdffszbaA3SYheitTkSjWVCR3PjaIEPPvgg2Oo9/fTT2b9/PwBffPEFOTk55OTkMHz4cI4ePUp+fj4TJkwgJyeHIUOGsHTpUgBeffVVhg4dypAhQ7j99tvrfZzvvvuOcePGMWzYMEaPHs3Ro0epqKjgl7/8JUOHDmX48OEsXrwYAK/Xy6233sqoUaPIzs7mX//6FwBTp06ltLSUMWPG8Pe//53bbruNBQsWkJOTQ3l5Ob179+bgwYMAvPDCC2RnZzNs2DCuuOIKAO655x4eeOABAH788UfOPPNMRo4cyfjx49m0aRMAs2bN4sYbb2TcuHH07duXN998E4A77riDpUuXkpOTw0MPPcT69esZPXo0OTk5ZGdns3Xr1tZ4eZokLXYhWtug8yGlB/TIbfRmf/pgPRv2Fof3obulcPd5g1t8v5NPPplvv/0WpRTPPPMM//jHP/jnP//JAw88wBNPPMFJJ51ESUkJcXFxzJ07lylTpnDXXXfh9XopKytj79693H777axcuZK0tDQmT57Mu+++y/nnnx98DJfLxSWXXMLrr7/OqFGjKC4uJj4+nkceeQSAtWvXsmnTJiZPnsyWLVt44YUXSE1N5bvvvqOyspKTTjqJyZMn8/7775OUlMTq1asB6Ny5MytWrODxxx+v8ZzWr1/PX/7yF7766isyMjIoLCys87xnz57N008/Tb9+/Vi2bBnXXnstn332GQD5+fl8+eWXbNq0ialTp3LhhRdy33338cADD/Dhhx8CcMMNN3DTTTcxY8YMXC4XXq+3xfs+HCSxC9HaEjpC/8lmR9EieXl5XHLJJeTn5+NyuYL11CeddBI333wzM2bMYPr06fTo0YNRo0Zx5ZVX4na7Of/888nJyeGzzz7jlFNOITB9yIwZM1iyZEmNxL5582a6du3KqFGjAEhJSQHgyy+/5IYbbgBgwIAB9OrViy1btrBw4UJ++OGHYGv5yJEjbN26tdm13p999hkXXnghGRkZAHTsWPMXVElJCV9//TUXXVR1nGVlZWXw8vnnn4/FYmHQoEHBXzC1nXjiifzlL38hLy+P6dOn069fv2bFFm6S2IVobbuWwbzJcMG/YNilDd7sWFrWreWGG27g5ptvZurUqXz++efcc889gNH1cM4557BgwQLGjh3LokWLmDBhAkuWLGH+/PlcccUV3HrrrcEk3Ritdb1lfQ1NJa615rHHHmPKlCnH9JwaerwAn89Hhw4dgi3/2pxOZ5MxXn755YwZM4b58+czZcoUnnnmGU477bRjijcU0scuRGvbvcxY7ltrbhwtcOTIEbp37w4QHIQEow966NCh3H777eTm5rJp0yZ++uknOnXqxDXXXMNVV13FqlWrGDNmDF988QUHDx7E6/Xy6quvMnHixBqPMWDAAPbu3ct3330HwNGjR/F4PEyYMIGXX34ZgC1btrBr1y5OOOEEpkyZwlNPPYXb7Q6uKy0tbfZzmjRpEm+88QaHDhmD2LW7YlJSUujTpw///e9/ASN5r1mzptFtJicnc/To0eD/27dvp2/fvtx4441MnTqVH374odnxhZO02IVobVF+AFJZWRk9evQI/n/zzTdzzz33cNFFF9G9e3fGjh3Ljh07AHj44YdZvHgxVquVQYMGcdZZZ/Haa69x//33Y7fbSUpK4oUXXqBr16787W9/49RTT0Vrzdlnn820adNqPK7D4eD111/nhhtuoLy8nPj4eBYtWsS1117LnDlzGDp0KDabjeeeew6n08nVV1/Nzp07GTFiBFprMjMzeffdd5v9PAcPHsxdd93FxIkTsVqtDB8+nOeee67GbV5++WV+/etfc++99+J2u7n00ksZNmxYg9vMzs7GZrMxbNgwZs2aRUVFBS+99BJ2u50uXbrwxz/+sdnxhVOLzqAULrm5uVpOtCFixtePwcLfw9jr4My/1li1ceNGBg4caFJgIprV995QSq3UWjc+Co90xQghRLsjiV2I1mZ1GMu0XubGIWKG9LEL0doGT4fOQ6D7CLMjETFCErsQrS0p0/gTIkLC0hWjlPqtUmq9UmqdUupVpVRcOLYrRLuQvwbuSYWvHjU7EhEjQk7sSqnuwI1ArtZ6CGAFGj4KQ4hYs2OJsSyp/2hFIcItXIOnNiBeKWUDEoC9YdquECJKNDX9LcAf//hHFi1adEzb//zzzzn33HOP6b7NsWnTpuDkZT/++COPPvooAwcOZMaMGbz//vvcd999jd5/3Lhxx/zYzz33HHv3Ri4thtzHrrXeo5R6AONk1uXAQq31wtq3U0rNBmYDZGVlhfqwQogI0VqjtWbBggVN3vbPf/5zBCI6Nu+++y7Tpk3jT3/6EwBPPvkkH330UXCumalTpzZ6/6+//vqYH/u5555jyJAhdOvW7Zi30RLh6IpJA6YBfYBuQKJS6ue1b6e1nqu1ztVa5wYmBhIiJphwEGBLPPjggwwZMoQhQ4bw8MMPA8Y84wMHDuTaa69lxIgR7N69u8b0t//3f//HgAEDOOOMM7jsssuC097OmjUrOElX7969ufvuuxkxYgRDhw4NToG7fPlyxo0bx/Dhwxk3bhybN29uND6v18stt9zC0KFDyc7O5rHHjBOWfPrppwwfPpyhQ4dy5ZVXBifsWrlyJRMnTmTkyJFMmTKF/Px8FixYwMMPP8wzzzzDqaeeypw5c9i+fTtTp07loYce4rnnnuP6668HYP/+/VxwwQUMGzaMYcOGBRN69ROS3H///cHpg+++++4a++yaa65h8ODBTJ48mfLyct58801WrFjBjBkzglMJ33HHHQwaNIjs7GxuueWW0F/EWsJRFXM6sENrXQCglHobGAe8FIZtC9H2+TzGssvQpm/77Dl1rxt8Poy+Blxl8PJFddfnXA7DZ0DpIXjjFzXX/XJ+ow+3cuVKnn32WZYtW4bWmjFjxjBx4kTS0tLYvHkzzz77LE8++WSN+6xYsYK33nqL77//Ho/Hw4gRIxg5cmS928/IyGDVqlU8+eSTPPDAAzzzzDMMGDCAJUuWYLPZWLRoEb/73e946623Goxx7ty57Nixg++//x6bzUZhYSEVFRXMmjWLTz/9lP79+/OLX/yCp556iuuuu44bbriB9957j8zMTF5//XXuuusu5s2bx5w5c0hKSgom0o8//pjFixeTkZFRY2qBG2+8kYkTJ/LOO+/g9XopKSmpEc/ChQvZunUry5cvR2vN1KlTWbJkCVlZWWzdupVXX32Vf//731x88cW89dZb/PznP+fxxx/ngQceIDc3l8LCQt555x02bdqEUqrJ7q1jEY7EvgsYq5RKwOiKmQTIfAFCBAy/Ao6fBJ2bkdgj7Msvv+SCCy4gMdE4Gcj06dNZunQpU6dOpVevXowdO7be+0ybNo34+HgAzjvvvAa3P336dABGjhzJ22+/DRgTjM2cOZOtW7eilApO6tWQRYsWMWfOHGw2I1117NiRNWvW0KdPH/r37w/AzJkzeeKJJzj99NNZt24dZ5xxBmC09rt27dqSXcJnn30WPFWe1WolNTW1xvqFCxeycOFChg8fDhjT/W7dupWsrCz69OlDTk5O8Dnv3LmzzvZTUlKIi4vj6quv5pxzzmmVcYVw9LEvU0q9CawCPMD3wNxQtytEu9GSOvbGWtiOhMbXJ6Y32UKvrbG5ogLJviX3qS0w1a3VasXjMX65/OEPf+DUU0/lnXfeYefOnZxyyilNxlh7ut3GpvYdPHgw33zzTbNjbCmtNXfeeSe/+tWvaly/c+fOGlP7Wq1WysvL69zfZrOxfPlyPv30U1577TUef/zx4Mk8wiUsVTFa67u11gO01kO01ldorSubvpcQMSJvhVHHvvD3ZkdSx4QJE3j33XcpKyujtLSUd955p8lzhJ588sl88MEHVFRUUFJSwvz5LfsyqT4lcO3ZFeszefJknn766eAXQ2FhIQMGDGDnzp1s27YNgBdffJGJEydywgknUFBQEEzsbreb9evXtyi+SZMm8dRTTwFGi7+4uOZZraZMmcK8efOCXTR79uzhwIEDjW6z+vS+JSUlHDlyhLPPPpuHH364wfnfQyFzxQjR2jb5E58n+to7I0aMYNasWYwePZoxY8Zw9dVXB7sYGjJq1CimTp3KsGHDmD59Orm5uXW6Kxpz2223ceedd3LSSSc169RxV199NVlZWcFzlb7yyivExcXx7LPPctFFFzF06FAsFgtz5szB4XDw5ptvcvvttzNs2DBycnJaXM3yyCOPsHjxYoYOHcrIkSPrfDFMnjyZyy+/nBNPPJGhQ4dy4YUX1piTvT6zZs1izpw55OTkcPToUc4991yys7OZOHEiDz30UIviaw6ZtleI1vbJH+GrRyD3Kjj3wRqr2uq0vSUlJSQlJVFWVsaECROYO3cuI0bIXDjhFMq0vTJXjBCtzedvlQaqY9qB2bNns2HDBioqKpg5c6Yk9SgjiV2I1hZI7NqcM9a3hldeecXsEEQjpI9diNYWaKn3Orne1WZ0h4roFup7QlrsQrS2U+6EE6+Djn3qrIqLi+PQoUOkp6fXKekTsUlrzaFDh4iLO/ZJciWxC9HaEtMhPs3okrFYa6zq0aMHeXl5FBQUmBSciEZxcXE1TjDeUpLYhWhtWxbCKxfBoGlw8Qs1Vtnt9uAkVEKEi/SxC9Ha1huH0gcHUYVoZZLYhWhtgcFT7TM3DhEzJLEL0draYR27iG6S2IVobYH6demKEREiiV2I1hZI6IMaP0OPEOEiVTFCtLbznwTvo0bZoxARIC12IVpbXCrY46Gy8RkAhQiXsCR2pVQHpdSbSqlNSqmNSqkTw7FdIdqF71+Gv3aFF6ebHYmIEeHqinkE+FhrfaFSygEkhGm7QrR9P7xuLNvRJGAiuoWc2JVSKcAEYBaA1toFuELdrhDthpQ7iggLR1dMX6AAeFYp9b1S6hmlVP0nSxQiFgXLHeUAJREZ4UjsNmAE8JTWejhQCtxR+0ZKqdlKqRVKqRUy4ZGIKcEjT6UrRkRGOBJ7HpCntV7m//9NjERfg9Z6rtY6V2udm5nZzDO2C9EeBLpiRl1tbhwiZoTcx6613qeU2q2UOkFrvRmYBGwIPTQh2okrPzbmibHHmx2JiBHhqoq5AXjZXxGzHfhlmLYrRNtnc0L5YSgvgpSuZkcjYkBYErvWejXQ5JmzhYhJS/8Jn/4ZUnvCb9eZHY2IAXLkqQid1rBsLpQVmh1JdFrzmrGUScBEhEhiF6Hb9wN8dCu8M8fsSKJToCpG6thFhEhiF6GzOo1l75PMjSNaBVrqUu4oIkQSuwif1GM/+W67JkeeigiTaXtF6DwVxvLgVnPjiFZKQce+cOJ1ZkciYoQkdhE6Z7KxzP/B3DiilVTCiAiTrhgRuo59AQWdB5sdSfQqPQgHNpkdhYgR0mIXoVMKrA7wyqSe9Zp/C2xbBEW74G4pCRWtTxK7CN2Pi8FbCbu+MTuS6LTmVXCVGJe1Nr4IhWhF0hUjQuepNJZyAE79qu8XLVP3itYniV2ELtAFc0qd2ZoF1Kxfly8/EQGS2EXoAom9Q5a5cUQrnwcs9qrLQrQy6WMXoQskqx1LIPMEc2OJNlob5aDHTYLjTgWr3eyIRAyQxC5Cl9bHWG54D0ZfY24s0UYpuGOX2VGIGCNdMSJ0WWOg89CqA5VEXSUFkLcCvG6zIxExQFrsInQ+H1htUsdeH3cFvDMbXKVGLfstWyGpk9lRiXYubC12pZRVKfW9UurDcG1TtBHfPAZ7vzfOECRq8rqMLqrDPxn/y+CpiIBwdsXcBGwM4/ZEWyHdCw0LlDpaHcZSyh1FBIQlsSulegDnAM+EY3uijQkk9kteMjeOaBRI5LZAYpcWu2h94WqxPwzcBshhdbHI6wKLLbZP1Pzj4von+Qok9sDJSOTIUxEBISd2pdS5wAGt9combjdbKbVCKbWioKAg1IcV0cTrMlqiK541OxLzvHg+LLil/nUpPaDfGXDhszJwKiIiHC32k4CpSqmdwGvAaUqpOr/JtdZztda5WuvczMzMMDysiBq9/KfE++Zxc+MwU8YJkNCx7vXJneHm9TDhFhgyXUpCRUSEnNi11ndqrXtorXsDlwKfaa1/HnJkou0YcDYMuyy2yx0PbYOi3Q2vLz1odNdUFEcuJhGz5AAlETpXmTGAGsvVMdoLe1fVvf7IHnhxOix90OiuKdwe+dhEzAlrYtdaf661Pjec2xRtwIe/hXVvxm6L3dfIgKirBH78FMoO+W8r5Y6i9UmLXYQukNC9MVrK5/XPR59Yz8Bo7XJHLYldtD5J7CJ0PjekZsFNq82OxBzucmM5/ua66wJ164FyR6ljFxEgiV2EzuuGhLT6q0JiQeAMUvWRI0+FCSSxi9B5XZC/Bhb9qfH+5vYqLgWUBT6u5wxSVodRCtnrRJjxJnQeHPn4RMyRxC5Cl30JpPWGLx80umVijSMRRs6ChIy66zoPhuuXw8DzjIOUYvVXjYgoSewidMMuhdyrjMuxWPJYUWzMbll2sOHblBfBxg/h6L7IxSViliR2Ebqj+435xiE2Sx4PbDQSe332rob/TIGN78PrM2BPPbXuQoSZJHYRupd+Bl/cZ1yOxRa7t9rgae0xhvLDsPtb44QbIOWOIiIksYvQVW+lx2KLPVAV02di3dkb60zbK4ldtD5J7CJ0XhcMvgB+XwCpPcyOJvI8/tb4lL8YpwisTteatlfq2EUESGIXofN5wJ5gtEqVMjuayAu02CtL6rbI6xx5GoPloCLiJLGL0HldcHArLLgVjuSZHU3kZY2F/mfBs2fWff7OZOg2wqhl/+VH0PdUc2IUMUUSuwjdKXdA95GwfC6U7Dc7mshL7QFDLzQu1z4Ktc94mL0YugyBXuMgSc5FIFqfJHYRutwrjYNvIDarYgq3w09fGZcD/e21ucpgzWtwcFvk4hIxSxJ7G7d0awFfbWvkwJhIKNgMlUeNy7GY2Ne/AyvmGZdrt9g3fwxPnQQFG+GdX8GOz5u92Y/X7eO0Bz7H7ZV+edEytqZvIqLZQ59swaIUJx1fz+HskaA1PDEG+kww/o/lckeo22IvL4T966r+b8FcOhvyi9l+sJQj5W4ykpwhBiliSThOZt1TKbVYKbVRKbVeKXVTOAITzVNS6aG4wsRWss8LaHAkgcUem3XagWQ+6Y/QIavmuhCm7S2t9NRYCtFc4Wixe4D/p7VepZRKBlYqpT7RWm8Iw7ZFE0oqPHi1Ni+AQAu952i47BXz4jCTpxLiUmH8/6u7zldr2t4WHHkaSOglkthFC4XjZNb5WutV/stHgY1A91C3K5qnpNJDcbmJH/xAYrfazYvBbJ4KUFY4vLPuyaoDLfRjOPK0JNhij8FfQSIkYR08VUr1BoYDy8K5XVHTT4dK2V9cgdaakkoP5W4vLo9JA2yBwdKKYnj3Wvjpa3PiMNPY64yjTh8ZBhs/qLkupZsx/hCXCr9aCjkzmr3ZkmCLPQYHpEVIwjZ4qpRKAt4CfqO1Lq5n/WxgNkBWVlbt1aIFbnj1e3p2TOCBC4fh8/fCHK1wk27GAJsjEc57FFK7GxOBZY016rVjSWZ/42QbUHNCMIAB5xh/APFpLdpsVVeMtNhFy4Slxa6UsmMk9Ze11m/Xdxut9Vytda7WOjczUw7SCMWB4koKiis5Wq0ld6TcpFadIwFGzoQu2cb/sVgVs/NL2PWtcbmx0+R99wzkrWj2ZgMJXQZPRUuFoypGAf8BNmqtHww9JNGU4go3xRXuGn2vxRUmffjd5ZC3Elwlxv+xWMe+9J+w5H7jcu1yx+/+Aw9nG/tpwa2w+aNmb1aqYsSxCkeL/STgCuA0pdRq/9/ZYdiuqIfb66PM5aW43E1JtWRebFaLvXAHPHMa/PSN8X8sJnZPJThTqi5XV34Yin4yBlcttgarYio9Xsb+9VPm/5AfvC7Qx37UrC9t0WaF3Meutf4SiMEp/cwRSOBHyt01yuBM64oJdL04EiGug5G8Yo2nwug/P+t+Y86c6gJVMBarkdwbqGM/VOJiX3EFG/KPcE52V6B6VYwkdtEyMfgpbNsCXS6lLm+NZG7aQUqBROVIhDt+MicGs3kqwRYHY2bXXRdooSuL8aXXwJGnh8tc/qXxOrq9vmClU6lLErtoGZkrpo2p3uWSf6S82vUmffiljt1osducxpw5RbtqrvN5jJa6UmCxNNhiL/In9CJ/gq/eSpeqGNFS0mJvY6q30vcWVUvsZrXYA4ndYoe3roG+p8Dw5tdqtwsXv2C02J+fajz/85+oWpdxQlW546+WVPXF1xJosQcSfPV+9RIzp4wQbZK02NuY6gl8b5FRgZHgsJrXx95pEFw4DzJPgC3/g31rzYnDTJ0HQ/pxRqu9dlXMsEvgkheNy2m9IaFjvZsIdMEEltW7X+TI02O3p6ic+z7ahNdn4rQbJpDE3sZUT+B7isqxKOicEmdeVUxSJxjyM0jMMLpjYrGOfdWLRsmnLa7h+djBqGPfsrDeVUWlgRZ7za4Yh9Uic8WE4LNNB3j6ix/ZeajU7FAiShJ7G1O9L31vUTmJDhsp8Xbz6tiP7oftn4Or1JjoKhYT+/ybYeP7/hZ7rXLHRffAoyOMy0sfgg3v1ruJqha7sf8C/eqdU50yeBqCwAD0geJGDhxrhySxtzHVu2IKSipJirORGm83r8W+Ywm8MA2K88Fqi706dp/P+DKzxdXfFVNRDBVFxmVLw+WOgZZ6hdtHhdsbbLF3To6TcscQVHqML8gDRxv5JdUOyeBpG3Ok3E2Cw0qZy4vWkOi0kRJnI6+wzJyAglUxNujQC+I7mBOHWQJzw9icMOEWo6yxOu01qmLAn9jr7y8PtNQDlwMHn3VOiWPtniNhDztWBFrs+4slsYsoVlzupktKHLsKy/D4NEnOQFeMyVUxVgfM+tCcGMwUaKHb4uD40+uu93mqDtpq5MjTQFcMwOHSqoPPOqfEUenx4fH6sFnlB3ZLVSV26YoRUay4wkNKvJ2UeKNuPMlpIyXOTnG5B23GCTeCZwhyRP6xo4GnWov94La6k3z5fEZLHRo98rSozEWnZGfwcrArJsW4Tipjjk2wj/1obCV2abHXw+XxUeHxkhIXfQfdHCl3kxpvJyXORmGpiySn0cfu8vqocPuId1gjG1D1A5Q+vtNYnvHnyMZgpoQMuH6lUcb48R2wexnctKZqfc/RVdP1/nJBg1MuHC5zM6BLMgeOVnK4zE2Jy4PDZqFDgvEeLHF5SE2IvvdjtKuUrhgR8OAnW/jwh70sve1UjMkro8fRcjdZHRNI9bfYE502UuKNl7G4wh35xD7gHEg/HuyJsHd1Veu0jTpUUsnRCg+9MxKbdwerDTKOZ29RORnYcdSuisn9ZdXlBmrYvT5NcYWbvpmJLNtRGOxjT3LaSHQar22JTAR2TKqqYmIrsUtXTD2+3X6IvMPl7D0SfW+GI+VuUuJs1bpirMFfFqZUxqT1hv5TjATXDurY735/PZf9+9vmd2sd3Y/+6lFueupdlu8ubXw+9hXz4PuX61x9pNyN1tA73fgyCXTFJDqtVYldKmOOictb1RVjSlelSSSx1+Lx+tiYb5wAam1edFUjaG207FLi7cFknlQtyZsygHpgY9Uc41ZHmy531FqzfEch+UcqyDtc3vQdAIp2oT75A/FHd7Cr2Fc3sf93FvxronF59auw9r91NhGoiOmcEkeCw2p0xVR6SXLaSfYn9vZa8vifL3eweNOBVtt+oMVe5vLG1JejJPZathWUBPvl1u+NrsRe4fbh9mqjj71aV0z3DnEAbNlfEvmg1rwGb8w0LrfxxL6nqDw4yLZq1+Hm3cljfAFUageHKhW6dh27p7Lm1L31DJ4Gatg7JNhJS3Bw2N9iT6rWYm+Pif1ohZu/LdjIQ4u2tNpjVFY7F3BIlTGHd8L+DeAyqay4hSSx1xJopSfH2aKufjgwnUBKnD3Yr57stHFcZhLdO8Tz6cb9kQ/K666qiEnrDWm9Ih9DmHy/q6jey43yt9ArsfOedxwrxz4G1X/y+7xV4w4WG+i60/YeLjVe17QEBx0S7BRs5W64AAAgAElEQVSVuSl1eUh02khqx10xX249iMen+SHvSKsdQBQ4QAlCPEhp+b/hqRPh773aROMlXOc8PVMptVkptU0pdUc4tmmWdXuOkOiwcsbAzqzbcySq+uUCXS2p8fYag6dKKU4f2Ikvtx2kwh3hsjivy+hfBzjzr3Bp3T7ktmLVrsPE2S3k9krj+2a32I1k0bdrOrusWXzsHm5M0Rvg81Qrd6x/2t5AV0xagiPYYi+pMBJ7e+5jX7z5ADaLsa++2FzQKo/h8vjI8J/kPaRpBbZ9aiy9Lsj7LgyRta5wnPPUCjwBnAUMAi5TSg0KdbtmWbe3mMHdUsnukcrBEldUHdgQGBxNibcF+9gDH/xJAztT4fbx9Y8HIxuU19VuathX7Soiu0cHcnt3ZP3e4mZ9SboqjK6YQT0zOa1LJZatH4O7WsuwGUeeBqbq7ZBoD7bYSyo9JDttJDqN+0ZbV4zWmvdW7+GFb3biO4aZE30+zeLNBUwZ3IXOKU4Wb26dfnaX10fPjvFACCWPxXuhYCOMv8V4LQNJPoqFo9xxNLBNa70dQCn1GjAN2BCGbdfw1baDeFa9hNNds4ukJK4bOztNAuCEPW9j99acya04PotdmcYA1sDdr2PVVZUb+1OHUZCaTWrpTnoeWsqo/AKGZ3VgUFEKT9GVxz7byvC4ffQo/LpOPNu6nEOFoyPpRzfR9XDdb/HN3c7HbUsm88g6Oh/5vs76jd0vwmuNo3PR92QWr6uzfn2Py9AWG10Ll5NesplDJZVcZS2i79YNpJV7gGxjcG3bIsYd2MC1ji3sXfA5365MwmtxsrHHxQBkFXxBSnnNE0C4rYls7j4dgN4HFpFUUXWuzX6dkknP7MK6zHNIcFjpu+9j8vN28lNhKTaLhZyeqVhSuvJm5Vhy8g/T1aN4fel2hv70PFkHl/BT5inBbRUl9CYvYzwAg3e/gqp15GVhUj/2dhwL2seQ3XVb+weTB7EvbSRWbyUD97yBT9nZ2vUc3LbkGrer+RwVOzpNojSua3B9d72fM20r8Po0a/KO4PZXS+zpdibnjR+FLtjMmH2vMKZPOj0q4nGr3Sx9/hsScy9lXM5gCratZPuy+WhqJrFv40/llYonuK//YLodeZGzd96He3VHVGUx8ytz6GQbhcXiYe3S7Viy/gYo+r7xKHHuql8ECQUlnGeNJ9l5NmkJDk4uXkCCLmPk4TScy5fwK/smXOs3817HS5iW092YJdJdQV5ROXmHy/z7uQ95GSfjtFu5zDcfK17W5B0JfjkVJvVnb8cxzd7PtQU+J3ZPKSfsfZv8ogrW+cegXluWRHy/iRxKGYjTXUS//A/q3H9Px7EcTupHfGUBx+3/mHKXl6llB5lm78ZZHct4eks/nlmaRmJFPn0OLKpz/52Zp1ES352Usl1kHfyizvrtnSaT3q0PUzKLKN/4P7ZWppE9+Re4PD4m2dcx1rGehJWf821eEgA+ZWdDz0sB6HlwKallO2tsz21NYHP3n1FS6WHQstuYDDy6fwgXpwwlfsUbLHKNrvF8atudPp4jib3rPJ+C5MHknHw2PTsm1LlPOIUjsXcHdlf7Pw8YU/tGSqnZwGyArKysY3qgj9ft44p1c+lv2VPj+i+82dzr7gbAV84n6a4O1Vg/3zuae92dAFjjfJRUVTUAssuXyQTXI8yz/4Ox1tWMtfifQR4MTriXl5ftosKyhH86nq4Tzz3rOrFJZ/EL6/+4wP58nfW3rMsiT3fiWuuHTLO/Xmf9r9f25zAp3GJ7j/Ns79VZ/4u1w3Bh527b25xj+x8A59iB5dDNYiM1/lV6ZSTC529jW/0yt1mAYuOvSCcyY81QAJ60v8JY6/Ia287TGVyxeiAAL9hfYKy12jzqW4HMgVxf1pUBXVJ42jWXrruXEUyTm6GoYw637U2gj5pIKrmsnr+R2db9/M6+iq5Fq4Kbesd7Eve6MwDY6HyEeFWzHPIlzyTu9aRiwcf2uAfr7IOnPedynyeBFEr5wb/+/fWH+NA7lkftj/OCdzLf+gay3vlbrKoq6a7fuIF7PVfQlUOU4WSkZQtnOR7ABlQ/K+ll6+OxpfVg8KFl/M76EuwCdsEf7EAeXLCzE88cfzyLPpnPZfv/WSe+31V2ozwuixF9OrFzRxrsBPv8mwBY4PoN//ONNm64bmPwPgsczzDIUnUawbHA4PhslPozA7omczLv0Mt6IPg+vNMKCw+MZPZrJzCmTzpdFv8Nyg7SA+jh38bb3pO5150OwOUJ96B8leRUi/NFz+nc60lpZD+fV2c/V/dP94U85rXTlUN8419/QeB4qSL489dlzPPC8SqPRc6697/NfQ1veD3kqG28619/qh3YANnAx54buXd+J8ZZ1vGKo5741tv43JfDZMt3zHU8VDe+dUl8pwv58sz99Pj8j2QDpSdOw+XxMdqzhJMtC4KfDYBiHc/lPwwD4DH7K4y1fltje/t0GlesNjoe3k4pZp8jiyc2OCnU2dxjf4FlX37CG14Pw9Q23qvn+b603s2HvhPrPJ+nPOeRfMKEVk/sKtQ+ZKXURcAUrfXV/v+vAEZrrW9o6D65ubl6xYoVDa1uUIXbi7u8uObgFBh9lw7/ASWVR+ve0WIFe0Kd9fbVz2Pf8BZlMz7A+fm9aGcy3nE3kmA3vu9cljgqfcrobqivPtmeYGy7ofWORCM2T2X99d3B9RX1D8g4koz+2mrrHTaF0+r/aR/nPxuPuxy8bjS65mnUnP5Wrbus/smnAutdpcFBvd+/u5Yt+0t4fc44sv/6NScfn8FLvxjCRU99SVqig8OlLjxeTXpKPOsOePjfbyZgqd6hV3v/W2xgj69/XfX1WoOrnqoeq92Yh0VrVFkBSY8NovKUP+IecD5JT4+g/KxH8PUYTcK8iVSc+RCe/mejyg6iEzLAmUzi44PJy5jAaZvO4+1rRrDncDm3vrmGV64ZS6/0RCY+9C1nDetJ92Qb//psHUtvO40O8XYqvV627i9h6tzVzDmlH88v3cplIzpx0+n9asZnT8DhsOO0Wflp82oOvXwV5MxgWcIpPPjFbr7+3ZnE2Wv1eLpK6ryH4x12bPHG61FafBjQJDqM96HH5+Oj9Qe54a3NLLhxPIPSFWjN+H8s5pQTMrntzAFgsaHt8Zzx4BeM6+Ggc0o8r6/YxSe/mYjTbmnRfq5/vcOYNkH7wFWK3aqIsxnvw0qvF5e2G+t9XuP9VpvNaWzD5zHer9R8L5drB57AlAvuekpN7fHGc/C6653zvkw7mPLo1/RItnLywTe4w/4aB6/bwvRnNzCmRzz3XTCQMletz0AzPh8WpUhURndjpVZG+aSrxHguVkcjzzfO2Ke1n4/VTnx84jHP+6OUWqm1zm3qduFosecBPav93wPYG4bt1hFntxJnT2viRvUf3Vfv+om/hYm/JRlgat3WmMP/B3agsSMRm7O+MSGut8eDPR4FJMfVsz4utfH7x1XNyNivZzfe27iFFXuNL6JytxcciRR6nHSKT+Fn2Rnc8fZaOFDCb0/vX/cw98b2f1OvTXxTr103sNhwekpwYnyY4ncvhTGz4K59xGuf8WFKrbYdrxuHMw4PNkqI56iGEhLo2DGd1A4JjOyTybLth+iSGkfPLl3okGa0ep3AkL5p5GTt4qkvfkRrCz8bN5Dk1PpPbQfQ7fhsTvf+H1fF92V3YRld01LI9M//UvN5NP4eTkypud4GZHQ0Ek9RuQucGXh9mrxyG2lp6SRXe76nDejE+6v30jHJw7DjepKRkVH3AZraz02ur/mcnP4/gx2o701YfX183U02sb7m/eu2dpOBX088jr99tInhVuPxKysr6ereRXb5YayOkSTHN/A5aurz4U+TTsBps9Z6/Y7t+ba2cFTFfAf0U0r1UUo5gEuB98Ow3chxldb9FRCjBnUzEtfbq4zurnJ/K6fC7SPObmVqTjeSnTasFsWlo3s2uJ1WoZRxztCKI8Y85wDr3jReO4u16oTa3z4N6942LnvdWG2O4HMp85+0IsHfGh7bN53tB0v5bmchY/um13nIS0dloTUM69mBgV0bTuoAdquFvhlJbDtwlG37S+jXKSkMT9oQmDPmSLWTXmsNHRNrDlxPGtCZUpeX3YXlTBrYOWyP3xbMHNeb7B6p9OlkJOrKygrGeZZzxe672/wR0S0VcmLXWnuA64H/ARuBN7TW60PdbsR8cjf8tRvMPcXsSKJCIHkt3GDUxFd4AondS5zdQoLDxm1nDeCmSf3onNJYS6WVXL0ITvs9VBZXXTf/Zlj816r/V/wHNvjHLLwubA4j+ZW5vJT5BxMT/HPqBJK526sZ27duS/XcYV0Z1DWFX0/s26zwju+cxMb8o2w/WMLxncOf2Iv8lVGH/KfSS0+q2Xo+6fgMnDbjYz1pQKewPX5bEGe38v71J9Nv4mWcU/lXSqwdUIHyUktsTaAWlknAtNYLgAXh2FbE2fzJKbGen6wxqEtKHGkJ9uD84BX+Fnu520u83UiGV4w18SCk9OOMpao22diKeTXnQg/05QL43Fjtxmtc7vJS7vJiUQST36BuKSQ7bZS4PIzuUzexJzhsLLhpfLPD69cpifk/5PsvJzdx6+brEG98OQVKIw+V+BN7rRZ7vMPKGYM6k3+kgm4dIt8FEA2syZms170p81nRPjdYaXBWzfYqtp5tfbpmG8uOzWuRtXdKKQZ2TeHrH43KogqPD621kdgjPXNkfTa8bwzCDp8BN2+CBwcY12f0r7pNYDIyreHM+6DjUFhymFKXh9JKLwkOW3DWTqtFccqATuQXldMhIfR6/P6dk6tdDl+LPc5uwWGzBI8+LvS32Dsm1Y35wYtz8MVw12JqeR4/t36C6+jxWLQXHxYsltg6yF4Se/8zjQ9/zgyzI4kag/yJPSPJSWmlB5fXh9bGT13T/fA6FO4wEntyl6rr04+vumzxJ3alYMyvcLi9wMeUubyUuz11vqDuvzA7bEMs1fvVj8sMX2JXStEh3s6RciOhHyo1qrDSE+sOzjpssZXEautwZCP32p/ls6Jp2PHiU7aYmzsl1p5vXRYrjP11VemgICfLqJIZ07cjFR4vFS6jFDIqEntg8PTbp+HlC2HyX4zra7TY/ZOReT2wfz1O9xEsKjB46g32rwfE2a1h+zXSKz0Rm0XRvUN88KjgcAkclQpVXTFpcvKNOuz+MZXS8gqe90xm/sj/mBxR5EliF3WcPaQrn/x2AoO6pqB11eRjdeqxzRCXYgycHtgA+T9Axz6QmgUZ1erLZ/wXrngHKorgqXGotW+S4LAZg6euqrGC1uCwWTihSzJDuoe/oZAaX5XYC0tddEiwy3lQ6+Gw+08nWF5BPukUp2ebHFHkSVeMqMNiUfTrnMzSrca8M4X+SapaMyE2mzPF6GOvKDKS/IBzjL/qHP4657JCY2m1E++wUu72UF5Piz3cnpmZi6MVEm5qvIM9RcbBLodKK+sMnAqD3Wkk9rKKcnLVdo7fvxuYY25QESZf96JBge6Jw9GU2ONSAW1MzORsoFW8+hVY+mC187E6SHRYKa30BqfDbU1dU+PrlCGGQ4cEO0f8r8WhEle9/esCnA5jv1RUVPAz61JyNj5gckSRJ4ldNCjQ9RI4EURU9LGPnAW37TAqXho6YnDrQiO5B6ZpsNqJ93fFlLdyV0xr6hBvD9axF5a66hycJAy2nrlMcj3EWvpjV96qA9diiHTFiAYFEmDgRBBRkdid/kqT9OMgMbP+21gdRms92GK3k+Dviqlv8LSt6JBgp8zlpdLjpbDUxag+ktjr5UjggKM7HSstWPGiY6yGHSSxi0Y4/Yk80GKPijr2wz8ZByRNvL3qYKXaApMvpXSFaU9CtxEkOPZRUmkk9nhH23zbpyZUHaRUWOYiQ1rs9Ssp4BrLB6w6Og473pg76hSkK0Y0IthiL4uiqpiyg/DVw3Bwa8O3CbTY49OMeve0XsTbrf4jTz1tt8Xun8Rq58HSeueJEX6lB7jR9yKdK37EFqMt9ij4pIpoFUjs0VUV4+9Xf/USWPVi/bex2I3pZSuOwK5lUFFMgsNKqctDmbttd8UAbD9onEimYysM0LYL/kTucrm4x/0Ldpw+1+SAIk8Su2hQXK2umKjoY69+IFl984YDnPV3uG077F0N8ybDvrUkOG0cLnWjddXMjm1NYL6YHw8Yz1u6YhoQaKF73ewjHdWxj7nxmEASu2hQdA6eVkvsDVXFBE4mHayKcZBgtwZPCN1WW+yBE5iv3WOckq41SirbBX8VjE15OdvyLR13ts35CUMhiV00qHa5Y1R0xdirTRXcUB375o/g3euqVcXYaiTzqBgEPgaBk5os21FIr/SEsM733q74W+w2vMy0LSRtfd3TVrZ3kthFg+IcVYOnFgV2qzI5Ir8ZbxnLhub32b8OVr9Udcoyq6NGJUxbbbEnO21Y/C/BzBN7Y7FEyesRbZI6c2eft3jbOx4bXlQM1rGHlNiVUvcrpTYppX5QSr2jlOrQ9L1EWxE4p2VgLvbAVLemS0iD406DpC71rw+Ut7mMQUasjhrJvK0mdotFkRpvJ9Fh5cLcHk3fIVZZrLjiM6nEgQ0vFknsLfYJMERrnQ1sAe4MPSQRLexWhdXfKoyK/vWAn76B4yZBZv/611v9g4rdR8JFz0NylxrdL/H2tjl4CjChfybXnno8KXGxl6yazevm3EPzGK02Gi12W+ztq5De4VrrhdX+/Ra4MLRwRDRRShHvH3SMqsS+daFx5vdx19e/PtBCS+oMXYYAkOioqqBJdEbRc2mhRy4dbnYI0U/7OHX/c3xnuVha7GFwJfBRGLcnokBgADWqBhx3fAF5y8HnrX+9PcGomCncDts+BY+rXXTFiGYKDp76+Ln7LtR5D5scUOQ1mdiVUouUUuvq+ZtW7TZ3AR7g5Ua2M1sptUIptaKgoCA80YtWF2ipR8VRp7VZGkjQw2fAHbtgz0p4aTq4y2p2xbTROnbRTBYrGoVNeSi2pqFi8HzGTb7DtdanN7ZeKTUTOBeYpHXDJxjTWs8F5gLk5ubG7gkZ25hAYo+KUseAGW/Brq+bvl21aXsTHK7g1QnR9FxEq/ApG3a8XGVbAJsVnHCm2SFFVEhNF6XUmcDtwEStdVl4QhLRJD7YYo+iZNjvdOOvIXtWwlePgtN/YmmrnQRHVbdNVHUriVbhs9iw4uNK3octlphL7KH+vn4cSAY+UUqtVko9HYaYRBSJysTelJIDsOFdYwlgsQW7XywKnDF+sudY8L+zlvIPzyXY8FZNMRBDQq2KOb7pW4m2zBkYPG1Lid1arY7d6gClgt0vCQ5b9NTji1Zjj0/Gg82f2GOvKib2vspEi8RH8+BpQwJ17NkXw4T/B0CCM5DY29AXlDhmAzY+zlSLxoYHrLGX5mLvGYsWicrB06YEWmgdehpHqAIOqwWrRUlijxFddrzFeMvxWPHFZFdMG2qGCTMEW+xtKSE6EiClO+xbC1uMY+iUvztGSh1jg7LYsCkvV2a+BhNuMzuciJPELhoV6IIJzBvTJnQdBjdvgIIt8OFvg1fHO6zSYo8VVjt2vHgdycYXfYyRxC4aFWipt8kSQa+rxhnqEySxxwxltRNHJVcc/TfsWGp2OBEniV00KtgV05ZKBI/ug5cuNOaUsVadZei4zCT6ZiSaGJiIFCOxuzir+E3Y+73Z4URcG/q0CjMEB0/bUkvX64Jtn0BFUY0W+zMzc7ln6mATAxOR4r7qc65132T8E4OTgMlIkmhUmzxAqVorvfqHWurXY4fTbsWO/2jjGKyKib1nLFokOHjaFhP7qGtg5CxTQxHmsHzzGL9xfGX8Iy12IWpqk3XsgQ9yWu/gfOwixmxbxNm2deBDjjwVora4ttoVk3EC/PQ1dB4UPEhJxBCrnfQuveDqbWZHYgoZPBWNGt6zA6cN6MQJXZLNDqX5bE64fjkU74FvnzI7GmEGix18HrBYjL8YE3vPWLRIp5Q45s0aRWp8G/w563XXHEgVscNihaKf4IObIP8Hs6OJOEnson16fiocWB+TA2cCcKZAZQmsfA6O5JkdTcRJYhftU+CgFGmxx6YLnoIr/2dcjsEvd0nson0K1C7HYEWE8PN5jGUM1rGHJbErpW5RSmmlVOydNVZEJ6sDep0Mp/7O7EiEGVY+D+/+2rgsib3llFI9gTOAXaGHI0SYWB3QIQtSu5sdiTDDvh/g8A6wJ8Rkd1w4WuwPAbcBOgzbEiI8ug2DLR/Drm/NjkSYwWIDZyrclQ9ZY8yOJuJCSuxKqanAHq31mjDFI0R4XPISVB6FLf8zOxJhBosNfG6zozBNk51PSqlFQJd6Vt0F/A6Y3JwHUkrNBmYDZGVltSBEIY6B1sYHOwYrIgTG6+4ug7euhtP/FHNdck222LXWp2uth9T+A7YDfYA1SqmdQA9glVKqvi8BtNZztda5WuvczMzMcD4HIep680pjKYk9NsV1MJZr/wuuUnNjMcExDxdrrdcCnQL/+5N7rtb6YBjiEiI0h/xzhMTgwJkATv4NJHWGd+eAVapihGhfJLHHrkAfewweyxC2rzKtde9wbUuIkCWkQ1IXyJlhdiTCDBveh/dvMC7HYHectNhF+2SLg6RMiEsxOxJhhsD8MFZnTLbYJbGL9imtF+xbC/s3mB2JMEOglf7b9ZCYbm4sJpDELtqnUVcbywOS2GNSYBqBGK1ll8Qu2ievy1jGYP+qoOp1f/kiY17+GCOJXbRPC/9gLKUqJjbFpxnL/etAtaHTOoaJJHbRPpUWGMsYHDgTwIBzYPwtoOTUeEK0P9IVE7t87pj9YpfELtqn3uONZa9x5sYhzLHrW/jqEfBWmh2JKSSxi/bJHm/0rUqLPTZVFJsdgakksYt2SoP2QnG+2YEIMwTmh/nlx+bGYRJJ7KJ96jzEWFYUmRuHMEewjt1jbhwmkcQu2qfAB1rKHWNTYND0vzPNjcMkkthF+/Tdf4xl5VFz4xDmiEs1lmWHzI3DJJLYRfvUfaSxTIi9eUIE0HkQDL4AMvqbHYkpYm8GehEbzvgTjJwFHXqaHYkwi9dd1dceY2LzWYv2z2qHzNhsrQng8E+w6UOzozBNyF0xSqkblFKblVLrlVL/CEdQQggRkhithgkIqcWulDoVmAZka60rlVKdmrqPEEK0usCBaVMfNzcOk4TaYv81cJ/WuhJAa30g9JCEECJEgXLHGG25h5rY+wPjlVLLlFJfKKVGNXRDpdRspdQKpdSKgoKCEB9WCCEaERg0/fA35sZhkia7YpRSi4Au9ay6y3//NGAsMAp4QynVV2uta99Yaz0XmAuQm5tbZ70QQoSNPc7sCEzVZGLXWp/e0Dql1K+Bt/2JfLlSygdkANIkF0KYx5kMnYfGbLlrqF0x7wKnASil+gMO4GCoQQkhRMh8HqljP0bzgHlKqXWAC5hZXzeMEEJElM8HBRuNvxgUUmLXWruAn4cpFiGECI8YPB1edbH5O0UI0f5Z7DDuerOjMEVsf60JIdovq92YLyYGSWIXQrRP7jL4JjaPPJWuGCFE+3TOg9Atx+woTCGJXQjRPo26yuwITCNdMUII0c5IYhdCiHZGErsQQrQzktiFEKKdkcQuhBDtjCR2IYRoZySxCyFEOyOJXQgh2hllxiy7SqkC4KdjvHsG0Tnne7TGBdEbm8TVMtEaF0RvbO0trl5a68ymbmRKYg+FUmqF1jrX7Dhqi9a4IHpjk7haJlrjguiNLVbjkq4YIYRoZySxCyFEO9MWE/tcswNoQLTGBdEbm8TVMtEaF0RvbDEZV5vrYxdCCNG4tthiF0II0Yg2ldiVUmcqpTYrpbYppe4wMY6eSqnFSqmNSqn1Sqmb/Nffo5Tao5Ra7f8724TYdiql1voff4X/uo5KqU+UUlv9y7QIx3RCtX2yWilVrJT6jVn7Syk1Tyl1QCm1rtp19e4jZXjU/577QSk1IsJx3a+U2uR/7HeUUh381/dWSpVX23dPRziuBl87pdSd/v21WSk1JcJxvV4tpp1KqdX+6yO5vxrKD5F7j2mt28QfYAV+BPoCDmANMMikWLoCI/yXk4EtwCDgHuAWk/fTTiCj1nX/AO7wX74D+LvJr+M+oJdZ+wuYAIwA1jW1j4CzgY8ABYwFlkU4rsmAzX/579Xi6l39dibsr3pfO//nYA3gBPr4P7PWSMVVa/0/gT+asL8ayg8Re4+1pRb7aGCb1nq71toFvAZMMyMQrXW+1nqV//JRYCPQ3YxYmmka8Lz/8vPA+SbGMgn4UWt9rAeohUxrvQQorHV1Q/toGvCCNnwLdFBKdY1UXFrrhVprj//fb4EerfHYLY2rEdOA17TWlVrrHcA2jM9uRONSSingYuDV1njsxjSSHyL2HmtLib07sLva/3lEQTJVSvUGhgPL/Fdd7/85NS/SXR5+GliolFqplJrtv66z1jofjDcd0MmEuAIupeaHzez9FdDQPoqm992VGC27gD5Kqe+VUl8opcabEE99r1207K/xwH6t9dZq10V8f9XKDxF7j7WlxK7quc7Ukh6lVBLwFvAbrXUx8BRwHJAD5GP8FIy0k7TWI4CzgOuUUhNMiKFeSikHMBX4r/+qaNhfTYmK951S6i7AA7zsvyofyNJaDwduBl5RSqVEMKSGXruo2F/AZdRsQER8f9WTHxq8aT3XhbTP2lJizwN6Vvu/B7DXpFhQStkxXrSXtdZvA2it92utvVprH/BvWuknaGO01nv9ywPAO/4Y9gd+2vmXByIdl99ZwCqt9X5/jKbvr2oa2kemv++UUjOBc4EZ2t8p6+/qOOS/vBKjL7t/pGJq5LWLhv1lA6YDrweui/T+qi8/EMH3WFtK7N8B/ZRSffwtv0uB980IxN9/9x9go9b6wWrXV+8XuwBYV/u+rRxXolIqOXAZY+BtHcZ+mum/2UzgvUjGVU2NVpTZ+6uWhvbR+8Av/JULY4EjgZ/TkVbpcUIAAAEUSURBVKCUOhO4HZiqtS6rdn2mUsrqv9wX6Adsj2BcDb127wOXKqWcSqk+/riWRyouv9OBTVrrvMAVkdxfDeUHIvkei8Qocbj+MEaPt2B8295lYhwnY/xU+gFY7f87G3gRWOu//n2ga4Tj6otRkbAGWB/YR0A68Cmw1b/saMI+SwAOAanVrjNlf2F8ueQDbozW0lUN7SOMn8lP+N9za4HcCMe1DaP/NfA+e9p/25/5X+M1wCrgvAjH1eBrB9zl31+bgbMiGZf/+ueAObVuG8n91VB+iNh7TI48FUKIdqYtdcUIIYRoBknsQgjRzkhiF0KIdkYSuxBCtDOS2IUQop2RxC6EEO2MJHYhhGhnJLELIUQ78/8BJOmqATbzV3YAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "#!/usr/bin/python\n",
    "# -*- coding: utf-8 -*-\n",
    "\n",
    "print(__doc__)\n",
    "\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import time\n",
    "\n",
    "from sklearn.linear_model import Lasso\n",
    "from sklearn.metrics import r2_score\n",
    "\n",
    "# 用于产生稀疏数据\n",
    "np.random.seed(int(time.time()))\n",
    "# 生成系数数据，样本为50个，参数为200维\n",
    "n_samples, n_features = 50, 200\n",
    "# 基于高斯函数生成数据\n",
    "X = np.random.randn(n_samples, n_features)\n",
    "# 每个变量对应的系数\n",
    "coef = 3 * np.random.randn(n_features)\n",
    "# 变量的下标\n",
    "inds = np.arange(n_features)\n",
    "# 变量下标随机排列\n",
    "np.random.shuffle(inds)\n",
    "# 仅仅保留10个变量的系数，其他系数全部设置为0\n",
    "# 生成稀疏参数\n",
    "coef[inds[10:]] = 0\n",
    "# 得到目标值，y\n",
    "y = np.dot(X, coef)\n",
    "# 为y添加噪声\n",
    "y += 0.01 * np.random.normal((n_samples,))\n",
    "\n",
    "# 将数据分为训练集和测试集\n",
    "n_samples = X.shape[0]\n",
    "X_train, y_train = X[:n_samples // 2], y[:n_samples // 2]\n",
    "X_test, y_test = X[n_samples // 2:], y[n_samples // 2:]\n",
    "\n",
    "# Lasso 回归的参数\n",
    "alpha = 0.1\n",
    "lasso = Lasso(max_iter=10000, alpha=alpha)\n",
    "\n",
    "# 基于训练数据，得到的模型的测试结果\n",
    "# 这里使用的是坐标轴下降算法（coordinate descent）\n",
    "y_pred_lasso = lasso.fit(X_train, y_train).predict(X_test)\n",
    "\n",
    "# 这里是R2可决系数（coefficient of determination）\n",
    "# 回归平方和（RSS）在总变差（TSS）中所占的比重称为可决系数\n",
    "# 可决系数可以作为综合度量回归模型对样本观测值拟合优度的度量指标。\n",
    "# 可决系数越大，说明在总变差中由模型作出了解释的部分占的比重越大，模型拟合优度越好。\n",
    "# 反之可决系数小，说明模型对样本观测值的拟合程度越差。\n",
    "# R2可决系数最好的效果是1。\n",
    "r2_score_lasso = r2_score(y_test, y_pred_lasso)\n",
    "\n",
    "print(\"测试集上的R2可决系数 : %f\" % r2_score_lasso)\n",
    "\n",
    "plt.plot(lasso.coef_, label='Lasso coefficients')\n",
    "plt.plot(coef, '--', label='original coefficients')\n",
    "plt.legend(loc='best')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.5.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
